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Abstract 

This article is intended to provide a pedagogical account of issues related to, and recent work on, grav- 
itational waves from coalescing compact binaries (composed of neutron stars and/or black holes). These 
waves are the most promising for kilometer-size interferometric detectors such as LIGO and VIRGO. 
Topics discussed include: interferometric detectors and their noise; coalescing compact binaries and 
their gravitational waveforms; the technique of matched filtering for signal detection and measurement; 
waveform calculations in post-Newtonian theory and in the black-hole perturbation approach; and the 
accuracy of the post-Newtonian expansion. 

Contents: 1. Introduction. 2. Interferometric detectors. 3. Detector noise. 4. More about detector 
noise. 5. Coalescing compact binaries. 6. Waveform according to the quadrupole formula. 7. Matched 
filtering. 8. The signal-to-noise ratio. 9. Signal detection. 10. Signal measurement. 11. Waveform 
calculations: post-Newtonian theory. 12. Waveform to second post-Newtonian order. 13. Waveform 
calculations: perturbation theory. 14. Luminosity from the perturbation approach. 15. Accuracy of the 
post-Newtonian expansion. 16. Conclusion. 

1 Introduction 

The existence of gravitational waves is an unambiguous prediction of the theory of general relativity Q . Yet, 
despite efforts originating in the early nineteen sixties, gravitational waves have not been detected directly. 
Nevertheless, observation of the binary pulsar PSR 1913+16 convinces us that gravitational waves do exist, 
and that they are correctly described by Einstein's theory g{. That gravitational waves have not yet been 
detected on Earth is simply due to their incredible weakness: typical waves would produce in a bulk of matter 
a strain AL/L, where L is the extension of the matter, of order 10~ 21 (lj. Needless to say, to measure this 
effect is a great challenge for experimentalists. 



2 Interferometric detectors 

There is reasonable hope that gravitational waves will be detected within the next ten years, thanks to a 
new generation of detectors which use interferometry to monitor the small displacements induced by the 
passage of a gravitational wave. Two groups are currently involved in building large-scale interferometers: 
the American LIGO team, and the French-Italian VIRGO team. 

The LIGO (Laser Interferometer Gravitational-wave Observatory) project || involves two detectors, to 
be built in Hanford, Washington, and in Livingston, Louisiana. Construction has begun at both sites. Each 
interferometer has an armlength of approximately 4 km. LIGO should be completed by the turn of the 
century. 
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The VIRGO (so named after the galaxy cluster) project Q involves a single interferometer, to be built 
near Pisa, Italy, with an armlength of approximately 3 km. VIRGO should also be completed by the turn 
of the century. 

The basic idea behind interferometric detectors is the following 

The interferometer is composed of two long (4 km for LIGO) vacuum pipes forming the letter L. A laser 
beam is split in two at the corner of the L, and is sent into each arm of the interferometer. Each beam then 
bounces off a mass which is suspended at each end of the L (a mirror has been coated onto each mass). The 
light is finally recombined at the beam splitter, and its intensity is measured by a photodiode. 

When no gravitational wave is present at the interferometer, the length of each arm is so adjusted that 
when measured by the photodiode, the light's intensity is precisely zero (the recombined beams are arranged 
to be precisely out of phase). However, when a gravitational wave passes through the interferometer, the 
armlengths are no longer constant, and the recombined beams no longer precisely out of phase. More 
precisely, during the first half of its cycle the gravitational wave increases the length of one arm, and 
decreases the length of the other. During the second half cycle, the first arm is now shorter, and the second 
arm longer. The light's intensity therefore oscillates with the gravitational- wave frequency. The intensity is 
a measure of AL/L = h, where L denotes the interferometer armlength, and h the gravitational- wave field. 



3 Detector noise 

Interferometers are subject to various sources of noise which limit the detector's sensitivity to gravitational 
waves. The relative importance of each source depends on the frequency at which the interferometer oscillates 
§• 

At low frequencies (/ < 10 Hz) the detector's sensitivity is limited by seismic noise, which is due to 
the Earth's seismic activity. At frequencies larger than 10 Hz the seismic noise can be eliminated with 
sophisticated isolation stacks; these fail at low frequencies. 

At high frequencies (/ > 100 Hz) the detector's sensitivity is limited by photon shot noise, which is due 
to statistical errors in the counting of photons by the photodiode. This source of noise can be reduced by 
increasing the laser power, or making use of "light recycling" S. 

At intermediate frequencies (/ between 10 Hz and 100 Hz) the noise is dominated by thermal noise, which 
is due to spurious motions of thermal origin. For example, the suspended masses are thermally excited and 
vibrate with their normal-mode frequencies; this evidently affects the recombined laser beam. 

Interferometers are therefore broad-band detectors, with good sensitivity in the range |s| 

10 Hz < / < 1 000 Hz. (1) 

The required sensitivity for full-scale interferometers is approximately h n ~ 10 -22 at peak sensitivity — a 
tall order. [The subscript n stands for "noise level"; we will define h n (f) precisely below. A plot of h n (f), 
appropriate for an interferometric detector with "advanced" sensitivity, is given in Fig. 1.] For comparison, 
we may mention that the Caltech 40 m prototype has already achieved h n ~ 10 -19 at peak sensitivity 
(/ = 450 Hz). It is not implausible that improved technology and a factor of 100 in armlength will permit 
to reach the desired goal. 



4 More about detector noise 

The detector noise can be measured when no gravitational wave is present at the interferometer, the typical 
situation. Then the detector output s(t) — AL(t)/L is given by noise alone: 

s(t)=n(t), (2) 

where n(t) represents the noise. The noise is a random process E): the function n(t) takes purely ran- 
dom values. Consequently, the noise can only be studied using statistical methods. In the presence of a 
gravitational wave, Eq. (0) must be replaced by s(t) = h(t) + n(t), where h(t) is the gravitational- wave field. 



2 



10 



hjf) 



10" 



Iff' 



Iff' 



10 



Iff' 




I I I I 

10 100 1000 10 000 

f (Hz) 

Figure 1: Noise level in an interferometric detector with advanced sensitivity. 

The statistical properties of the noise can be determined by careful measurement. For example, the time 
average 

1 f+ T 

»(*)= 1™ ipr / "(*)<** (3) 

I — >00 Zl J_ rp 

can be constructed. This mean value can then be subtracted from n(i) and, without loss of generality, we 
can put n = 0. Also from measurements, the noise's autocorrelation function C n (r) can be constructed: 



C n (T)=n(t)n(t + T); (4) 

C n (0) gives the mean squared deviation of the noise with respect to the mean value. 

In the following we will assume that the noise is stationary, in the sense that its statistical properties 
do not depend on time j5j. This means, in particular, that the autocorrelation function does not depend 
explicitly on the origin of time t, but only on the variable r, as was expressed in Eq. (^). We shall also 
assume that the noise satisfies the ergodic hypothesis, so that time averages can be replaced with ensemble 
averages ||. Here, the noise is imagined to be drawn from a representative ensemble, and the probability 
that it takes a particular realization n(t) is given by a specified probability distribution. The statistical 
properties of the noise then refer to this (infinite dimensional) distribution function. 

In full generality, the statistical properties of the noise can only be summarized by constructing all the 
higher moments n ■ ■ ■ n. If, however, the noise is assumed to be Gaussian, in the sense that its probability 
distribution function is an infinite dimensional Gaussian distribution Q, then the autocorrelation function 
contains all the information. 

Real detector noise is neither strictly stationary nor strictly Gaussian. However, on a timescale of 
hours, which is long compared with typical gravitational- wave bursts, the noise appears stationary to a 
good approximation ||. And non-Gaussian components to the noise can be removed, to a large extent, by 
cross-correlating detector outputs from two widely separated interferometers || . It is therefore a satisfactory 
approximation to take the noise to be stationary and Gaussian. 

Under these assumptions the statistical properties of the detector noise are fully summarized by C„(r). 
It is convenient to work instead in the frequency domain, and to define (^] the noise's spectral density S n (f) 
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S n (f) = 2 J ' C n (r)e 2 ^ T dr. (5) 

The spectral density is defined for / > only; as C„(r) is a real and even function, the negative frequencies 
only duplicate the information contained in the positive frequencies. 

As n{t) is dimensionless, the spectral density has dimensions of time. By multiplying S n (f) with the 
frequency and taking the square root (since the spectral density represents the mean squared noise), one 
obtains the noise level h n {f): 

K{f) = y r fSjf)- (6) 
This gives the equivalent gravitational-wave amplitude which would make the interferometer oscillate at just 
the noise level El; this quantity was introduced in the preceding section (see Fig. 1). 

5 Coalescing compact binaries 

Coalescing compact binaries, composed of neutron stars and/or black holes, are the most promising source 
of gravitational waves for interferometric detectors jj], (?J . 

Consider the binary pulsar PSR 1913+16 This system consists of two neutron stars, each of 1.4 M Q , 
in orbital motion around each other. Its present orbital period P is approximately 8 hours, corresponding 
to orbital separations of about 5 x 10 5 times the total mass. (Here and throughout we use units such that 
G = c = 1.) Its present eccentricity is approximately equal to 0.6. 

The binary's orbital period is observed to decay at a rate dP/dt = —2 x 1CP 12 corresponding precisely 
to a loss of energy and angular momentum to gravitational waves q]. In a timescale of approximately 
10 8 years the orbital period will have decreased to less than a tenth of a second, corresponding to orbital 
separations smaller than one hundred times the total mass. In this time, the eccentricity will have been 
reduced to extremely small values (by the radiation reaction), so that the orbits are practically circular. The 
gravitational waves produced then have a frequency larger than 10 Hz, and the frequency keeps increasing 
as the system evolves. During this late stage of orbital evolution, the gravitational waves sweep through the 
frequency bandwith of interferometric detectors, and thus become visible. By the time the orbital separation 
becomes as small as a few times the total mass, the neutron stars begin to merge. The gravitational waves 
produced during the final merger cannot be detected by interferometric detectors, at least in the broad-band 
configuration described above ||: the gravitational-wave frequency is then larger than 1 000 Hz, for which 
the detector noise is large. 

Of course, it would be foolish to wait 10 8 years in order for PSR 1913+16 to produce gravitational waves 
with appropriate frequency. Fortunately, interferometric detectors will be sensitive enough to monitor binary 
coalescences occurring in quite a large volume of the universe, approximately 10 7 Mpc 3 (corresponding to a 
radius of 200 Mpc ||). It has been estimated || that as many as 100 coalescences could occur every year in 
such a volume (this includes coalescences of black- hole systems as well) . This potentially large event rate is 
one of the factors that make gravitational waves from coalescing compact binaries especially promising. 

The other factor comes from the fact that compact binaries are extremely clean astrophysical systems. It 
can be estimated |lCj that tidal interactions between the two stars are completely negligible, up to the point 
where the objects are about to merge. In particular, no mass transfer occurs. The system can therefore be 
modeled, to extremely good accuracy, as that of two point masses with a limited number of internal properties 
(such as mass, spin, and quadrupole moment). The challenge in modeling coalescing compact binaries resides 
in formulating and solving the equations of motion and wave generation for a general relativistic two-body 
problem JTT) . 

6 Waveform according to the quadrupole formula 

At the crudest level, the gravitational-wave signal corresponding to an inspiraling binary system can be 
calculated by (i) assuming that the orbital motion is Newtonian (with the effects of radiation reaction 
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incorporated) and (ii) using the standard quadrupole formula for wave generation [ p"2[ . As motivated above, 
we may also assume that the orbits are circular. 

We define h(t) to be the gravitational-wave signal. This is given by a linear combination, appropriate for 
interferometric detectors, of the two fundamental polarizations, h+ and h x , of the gravitational- wave field. 
At this level of approximation, the signal is given by Q 

h(t) = Q(angles) (M/r) (TrMf) 2/3 cos$(t). (7) 

As expected, the signal decays as the inverse power of r, the distance to the source. 

In Eq. (0), Q is a function of all the angles relevant to the problem: position of the source in the sky, 
orientation of the orbital plane, position and orientation of the detector on Earth. The parameter M is 
called the chirp mass and represents a particular combination of the masses, given by 

M = (mim 2 ) 3/5 /(mi + m 2 ) 1/5 . (8) 

The waveform depends on the chirp mass only, and not on any other combination of the masses. The 
symbol / represents the gravitational-wave frequency, which is equal to twice the orbital frequency. Because 
the system loses energy and angular momentum to gravitational waves, the frequency is not constant, but 
increases in time according to 

I = ^' M »' w < 9 > 

As a consequence, Eq. (jfy shows that the amplitude of the signal, which is proportional to (nAif) 2 ^ 3 , also 
increases with time. A signal which increases both in frequency and in amplitude is known as a chirp, and 
this is the origin of the term "chirp mass". Finally, the phase function $(i) is given by 



$(*) = J 2-!rf(t')dt'. 



(10) 



Because the frequency is not a constant, the phase accumulates nonlinearly with time. 

For a system of two neutron stars, the gravitational- wave signal undergoes approximately 16 000 oscil- 
lations as it sweeps through the frequency bandwith of an interferometric detector. The timescale for the 
frequency sweep is approximately 15 minutes. The orbital separation ranges from approximately 180 to 10 
times the total mass M = mi + rri2, and the orbital velocity 

v = (ttM/) 1 / 3 (11) 

ranges from approximately 0.1 to 0.4. This indicates that relativistic corrections must be inserted in Eqs. (0) 
and (0) in order to obtain a satisfactory degree of accuracy. 



7 Matched filtering 

How does one go about finding a gravitational-wave signal in a noisy data stream, when typically the signal 
is not very strong? And once the signal is found, how does one go about extracting the information it 
contains? For signals of precisely known form, one goes about this using the technique of matched filtering 
. Signals from inspiraling compact binaries, since they can be calculated with high precision, belong to 
this class. The basic idea behind matched filtering is to use our knowledge about the signal in order to go 
find it in the data stream, after the noisy frequencies have been filtered out. 

Suppose a signal of known form h(t; jt) is present in the data stream. Here, the vector /2 collectively 
denotes all the parameters characterizing the signal. In the case of inspiraling binaries, these would be the 
time of arrival, the initial phase, the distance to the source, the position angles, the chirp mass, and other 
parameters to be introduced in Sec. 12. The detector output is given by 

s(t) = h(t;p) + n(t), (12) 
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where n(t) is the noise, whose statistical properties are fully summarized by the spectral density S n (f), as 
discussed in Sec. 4. 

The first step in matched filtering |l4| is to pass s(i) through a linear filter which removes the noisy 
frequencies. The idea here is to use our knowledge about the detector noise contained in S n (f) in order to 
discard that part of the detector output for which the detector noise is large. The output is imagined to be 
decomposed into Fourier modes according to s(t) = J s(f)e~ 2 '* l f t df; the filter suppresses the modes s(f) for 
which S n (f) is large. 

In mathematical terminology, a linear filter is a linear operation on a function s(i). This operation can 
always be written as 

s[t) -> / w(t-t')s(t') dt', (13) 



where wit — t') is the filter function. The filter which removes the noisy frequencies is the one such that 

Ml) - ssiwr (14) 

the factor of 2 is conventional. 

The next step in matched filtering consists of computing the overlap integral between the filtered output 
( |l3| ) and the known signal h(t; Atrial)- The true value of the source parameters is not known prior to 
the measurement. This operation must therefore be repeated for a large number of trial values /Itriai; the 
corresponding signals h(t; Atrial) are known as templates. The overlap integral defines the function 5(/2triai) 
given by 

S(PtTial) = f h(t; p t riai)w(t - t')s{t') dt' dt 

= (^Gutrial) | s). (15) 

To obtain the second line we have inserted the Fourier decompositions of h, w, and s, and carried out the 
integrations over time. The inner product (-|-) is defined as 



where a and b are arbitrary functions of time, with Fourier transforms d and b. 

A quantity analogous to S'(/ttriai) can be defined for the noise alone, by replacing s to the right of Eq. ( |l5| ) 
by n. Operationally, this amounts to filtering the detector output when a gravitational-wave signal is known 
not to be present. Because the noise is a random function, the integrals J h(t; fttiia.l)w(t — t')n(t') dt' dt are 
random also. And because the noise is assumed to have zero mean, to take an average over all possible 
realizations of the noise (by repeating the measurements many times) would yield a zero value. None of 
these quantities would be especially useful. We are therefore led to consider the root-mean-square average 
of these integrals, 



NQltnai) = rms / h(t; /itr>alV(* - t')n(t') dt' dt (17) 

= \J (M/Wl) I MAWl)) = p(Mtrial)- 

To go from the first to the second line requires some machinery which will not be presented in this review. 
We refer the reader to Ref. [O] for the missing steps. 
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8 The signal-to-noise ratio 



The signal-to-noise ratio is defined to be the ratio of S over TV: 

amr- ^ (Mmw) I s) 

SNR(/X tria l) = —; r — . (18) 

It is clear that SNR(/2 t riai) is a random function, since s(t) is itself a random function, being the superposition 
of signal plus noise. Its expectation value, or average over all possible realizations of the noise, is zero in the 
absence of signal (since n = 0), and 

cm) M 7 (ft (Atrial) | Hp)) na . 

SNR(^trial) = TH \ L (19) 

PitHrial) 

in the presence of the signal h(t; p). It is important to notice that in Eq. (pl|), the numerator is the overlap 
integral between the true signal h(t; p) and the templates h(t] Atrial)- In the absence or presence of a signal, 
the variance in the signal-to-noise ratio is precisely equal to unity fl4|| , independently of the values of p 
and Atrial- This fully summarizes the statistical properties of the signal-to- noise ratio, since SNR(/7 tr iai) 
is a Gaussian random function. [This can be seen from the fact that s(t) itself is Gaussian, being the 
superposition of signal plus Gaussian noise.] 

It is intuitively clear that choosing a template with /Stria] = U will produce the largest possible expectation 
value of the signal-to-noise ratio. This statement, known as Wiener's theorem |i3|] , can easily be shown to 
be true by applying Schwarz's inequality to the right-hand side of Eq. (|l9|). We have 

max{sNR(/Ztriai)} = SNR(^) = p{p) = (h(p) \ h(p)). (20) 

The fact that the signal-to-noise ratio has a variance of unity indicates that a signal can be concluded to be 
present only if p(jT) is significantly larger than 1. We will come back to this point in the next section. 

The true value of the source parameters can therefore be determined by maximizing SNR(/2 tr i a i) over all 
possible values of the trial parameters. The fact that the signal-to-noise ratio has a variance of unity implies 
that this determination can only have a limited degree of accuracy. The statistical errors decrease with 
increasing p(p); since this is proportional to the signal's amplitude, a stronger signal gives better accuracy. 

Maximizing the signal-to-noise ratio is essentially equivalent to maximizing the overlap integral 

(ft(/Ztrial) | Hp)). 

It is easy to see how the choice of the parameters Atrial can affect the overlap integral. Consider a toy 
waveform with three parameters: arrival time, initial phase, and chirp mass. A mismatch in the arrival 
times clearly reduces the overlap integral: the signal and the template, taken to be functions of time, might 
have support in entirely different regions of the time axis, leading to a vanishing overlap. Supposing that 
the arrival times are matched, a mismatch in the initial phases can also reduce the overlap integral: the 
signal and the template might be out of phase with each other, leading to an oscillating integrand and a 
vanishingly small overlap. Finally, supposing that both the arrival times and initial phases are matched, a 
mismatch in the chirp masses would also reduce the overlap integral. This is because the chirp mass governs 
the rate at which the signal's frequency changes with time; cf. Eq. (||). Signal and template, starting at the 
same time with the same phase, might thereafter go out of phase, thereby reducing the overlap. 



9 Signal detection 

The first order of business when analysing the output of a gravitational- wave detector is to decide whether 
or not a signal is present. Here we assume that the signal must be of a specific form, corresponding to 
a coalescing binary system. In this section we discuss signal detection — how the technique of matched 
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filtering can be employed to find the signal in the noisy data stream. In the next section we will discuss 
signal measurement — how matched filtering is used to estimate the value of the source parameters once the 
signal has been found. 

As mentioned in the previous section, a signal can be concluded to be present if the maximum value 
of the signal-to-noise ratio, SNR = max{SNR(/2 tr i a i)}, is significantly larger than unity. In fact, there 
exists a threshold value SNR* such that a signal is concluded to present, with a certain confidence level, if 
SNR > SNR*. To figure out how large this threshold must be is a standard application of the statistical 
theory of signal detection which was developed largely for the purpose of detecting radar signals. The 
theory can easily be taken over to the case gravitational- wave signals [[IB). To go into the detail of this 
theory would be outside the scope of this review. We shall simply state that SNR* is fixed by selecting a 
small, acceptable value for the probability that a signal would falsely be concluded to be present — the false 
alarm probability. (This is the Neyman-Pearson criterion || , which is more precisely formulated in terms of 
the likelihood ratio, the ratio of the probability that a signal is present to the probability that it is absent.) 
Once SNR* is fixed, the level of confidence that a signal is indeed present increases with SNR > SNR*. A 
typical ballpark value for the threshold is SNR* = 6. 

In the preceding paragraph it was assumed that the signal-to-noise ratio SNR is computed using template 
waveforms which are functions of the parameters pi introduced in Sec. 7. These parameters have direct 
physical meaning; they include the chirp mass and other meaningful parameters to be introduced in Sec. 12. 
However, since these parameters arc not estimated during the detection stage of the data analysis (they are 
estimated only after a signal has been found), there is no particular need to parametrize the templates with 
ft. In fact, it may be desirable, in order to minimize the computational effort, to parametrize the signal in a 
completely different way. The new parameters, a, would then have no particular physical significance. What 
is required is that the new templates, h(t;ct), reproduce the behaviour of the expected gravitational-wave 
signal. In other words, the templates h(t; ft) and h(t;a) should span the same "signal space", but h{t;ct) 
should do so most efficiently. Unlike h(t; /2), h(t; a) need not be derived from the field equations of general 
relativity. With these new detection templates, the signal-to-noise ratio is defined as 

SNR = max{SNR(d?)}, (21) 

where SNR(d-) is defined as in Eq. (]l|). As before, a signal is concluded to be present if SNR > SNR* . 

10 Signal measurement 

Once the detection templates h(t; a) have been used to conclude that a signal is present in the data stream, 
they are replaced with the measurement templates h(t; Atrial) in order to estimate the value of the physical 
parameters ft. 

The procedure to estimate the source parameters was explained in Sec. 8. As was mentioned, the idea is 
to maximize over all possible values of the trial parameters the overlap integral 

(h(tf) | MMtrial))- 

We already have discussed the effect of a mismatch in the value of the parameters. What remains to be 
discussed is the effect of a mismatch in the functional form of the template with respect to that of the true 
waveform. 

The true waveform is governed by the exact laws of general relativity. The template, on the other hand, 
is necessarily constructed using an approximation to the exact laws. (That approximations must be made in 
analytic calculations is obvious; in a numerical treatment the approximation resides in the finite differencing 
of the field equations.) This, clearly, must have an effect on our strategy for extracting the information 
contained in the gravitational-wave signal. This can be seen simply from the fact that Wiener's theorem, as 
stated in Sec. 8, strictly requires signal and template to have the same functional form; they are allowed to 
differ only in the value of their parameters. 
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Suppose that gravitational waves coming from a given source are received without noise, so that the 
true signal h(t; ft) is measured accurately. Suppose also that the true value /2 of the source parameters 
is known (God has spoken). Then a computation of the overlap integral with Atrial = M but with an 
approximate template h(t; Atrial) will not yield the maximum possible value for SNR(/2 tr iai)- This is because 
the approximation differs from the true signal, both in amplitude and in phase. Since both signal and 
template undergo a large number of oscillations (recall that for a system of two neutron stars, this number 
is approximately 16 000), the overlap integral is most sensitive to phase differences: a slight phase lag causes 
the integrand to oscillate, thereby severely reducing the signal-to-noise ratio with respect to its maximum 
possible value. 

A gravitational-wave astronomer doesn't know before the measurement the true value of the source pa- 
rameters, and must work with an approximation to the true general-relativistic waveform. We have seen 
that the phase lag occurring between the true signal and the approximate template when the parameters 
are matched reduces the signal-to-noise ratio from its maximum possible value. It follows that maximizing 
SNR(^Ztriai) with approximate templates introduce systematic errors into the estimation of the source param- 
eters: evaluating the signal-to-noise ratio with Atrial = will return, for some <5/2, a number larger than 
SNR(/2). The systematic errors are precisely the value of 6(1 for which the signal-to-noise ratio is largest. If 
the templates are a poor approximation to the true signal, then the systematic errors will be larger than the 
statistical errors arising because SNR(/2 tr iai) is a random function (see Sec. 8). 

We therefore appreciate the need for constructing measurement templates which are as accurate as pos- 
sible, especially in phase [TtJ . The requirement is that the systematic errors in the estimated parameters 
must be smaller than the statistical errors. An estimate for the required degree of accuracy comes from the 
observation that the overlap integral will be significantly reduced if the template loses phase by as much as 
one wave cycle with respect to the true signal. Since the total number of wave cycles is approximately 16 
000, we have 

accuracy ~ 167)00 ~ 10 " (22) 

Since the orbital velocity v is of order 10 _1 when the gravitational-wave frequency is in the relevant band- 
width, relativistic corrections at least of order v 4 are required to improve the quadrupole-formula expression 
given in Sec. 6. As we shall see, this is an underestimate. 

11 Waveform calculations: post-Newtonian theory 

We have seen that accurate measurement templates are required to make the most of the gravitational- 
wave signals originating from coalescing compact binaries. We also have seen that the quadrupole-formula 
waveform, Eq. (R), is not sufficiently accurate; relativistic corrections at least of order w 4 are required. How 
does one go about calculating these? 

A possible line of approach is to use a slow-motion approximation to the equations of general relativity. 
This is based on the requirement that if v is a typical velocity inside the matter source, then 

u<l. (23) 

We shall call this approximation "post-Newtonian theory" |ll|. We point out that for binary systems, 
post-Newtonian theory makes no assumption regarding the relative size of the two masses. This is to be 
contrasted with the perturbation approach, discussed in Sec. 13, in which the mass ratio is assumed to be 
small, but no restriction is put on v. 



In post-Newtonian theory 19 , one starts by defining fields h a/3 as 

h aP = v^99 aP - V af3 , (24) 

where g a P is the inverse of the true metric g a p with determinant g, and rf 1 ^ is the metric of Minkowski 
spacetime. When the harmonic gauge conditions 

d f3 h af3 = (25) 
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are imposed, the exact Einstein field equations reduce to 

Uh af3 = 16n{-g)T al3 + A a/3 . (26) 

Here, □ = rj a ^d a dp is the flat-spacetime wave operator, T a @ is the stress-energy tensor of the source, and 
K a<3 is nonlinear in h a P and represents an effective stress-energy tensor for the gravitational field. 

The detailed way in which one solves these equations is quite complicated, and will not be described 
here. The essential ideas are these [|l9): 

One first integrates the equations in the near zone (r < A, where r is the flat-space radius and A the 
gravitational wavelength) assuming slow motion, or dh a/3 /dt <C dh a P /dx l . One does this by iterations: the 
nonlinear terms in Eq. (|26| ) are first neglected, and the resulting linear equations integrated. These solutions 
are then used as input for the next iteration. This process is continued until the desired degree of accuracy 
is obtained. This is the standard post- Newtonian approach jnj. 

One next integrates the equations everywhere in the vacuum region outside the source. This is done 
once again by iterations, assuming h a/3 -C 1, but assuming nothing about the relative size of dh a @ /dt with 
respect to the spatial derivatives. This is because the vacuum region contains the wave zone, in which the 
field propagates with the speed of light; a slow-motion assumption would therefore not do for the field itself. 
This is the post-Minkowskian approach [[ll] . 

Using the post-Minkowskian approach one constructs, by successive approximations, the most general 
solution to the Einstein equations outside the source. This is characterized by two infinite sets of arbitrary 
multipole moments [ pp| , the mass moments Mg m (t — r) and the current moments Ji m (t — r), were t and m 
are the standard spherical-harmonic indices. (In practice, the fields h a P are equivalently expressed in terms 
of symmetric-trace-free moments, not spherical-harmonic moments.) The general solution is then matched 
to the near-zone solution in the region of common validity, and the multipole moments are thus determined. 

Finally, one expresses the radiation field — the time- varying, 0(l/r) part of the gravitational field — 
in terms of the derivatives of the mass and current multipole moments [ pp| . This gives the gravitational 
waveform. The gravitational-wave luminosity dE/dt can also be obtained from the radiation field. 



12 Waveform to second post-Newtonian order 

To date, the post-Newtonian calculation of the waveform has been carried out accurately through second 
post-Newtonian order — 0(v 4 ) — beyond the leading-order, quadrupolc-formula expressions given in Sec. 6. 
The complete waveform will not be displayed here. Instead, we will focus solely on the waveform's phasing. 

The phasing of the waves can be determined from df /dt, the rate of change of the gravitational- wave 
frequency. This can be expressed as 

df_ = dE/dt 

dt dE/df [ ' 

where dE/dt is (minus) the gravitational- wave luminosity, and dE/df relates orbital energy to orbital fre- 
quency (/ is twice the orbital frequency; the orbits are assumed to be circular). Both these quantities can 
be expanded in powers of 

v^{nMff/\ (28) 

with leading-order terms jl2) 

dE\ 32 9 , n f dE\ 7r , 

where the subscript N stands for "Newtonian". We have introduced the reduced mass fi, the total mass M, 
and the mass ratio 77 as 

M= ; , M = m 1 +m 2 , ^ = 77- (30) 

mi + rri2 M 
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Notice that 77 is restricted to the interval < 77 < 1/4, with 77 = 1/4 for mi = m.2- It is easy to check that 
Eqs. ( p7| ) and (|2^) reproduce Eq. (^|) above; the chirp mass can be expressed as M. — rf/ 5 M. 

For simplicity we will, in the following, focus on the quantity dE/dt. To second post- Newtonian order, 
the luminosity takes the form 



dE 
~dt 



dE\ 
dt ) N 



1 



( 



1247 35 



V 336 12 



+ — 77 \v 2 + (4tt - SO)tj 3 



/ 44711 9271 
V 9072 ~ "504 



V 



65 
18 



r/ 2 + SS ) v 4 + 



(31) 



Here, the terms SO and SS are due to spin-orbit and spin-spin interactions, respectively 21 ; these occur if 
the masses mi and 777,2 are rotating. Let Si and S2 be the spin angular momentum of each mass, and define 



the dimensionless c 
momentum. Then [21| 



is the spin-orbit term, and 



uantities \ a = S a /m a 2 , for a — {1,2}. Let also L be the direction of orbital angular 



SO = 



ll(m a /M) 2 + 577 



Xa 



V 



SS = -M 103x1 • X2 ~ 289L • xi L ■ f 2 



(32) 



(33) 



is the spin-spin term. 

In Eq. ([jl]), the leading-order term was first calculated in 1963 by Peters and Mathews 0] using the 
standard quadrupole formula. The first post-Newtonian correction, at order v 2 , was calculated in 1976 by 
Wagoner and Will |2j|. The 47™ 3 term is due to wave propagation effects: As the waves propagate out 
of the near zone they are scattered by the curvature of spacctime, and this modifies both the amplitude 
and the phase of the waveform; this "tail term" was first calculated in 1993 by this author (24], and then 
independently by Wiseman pj| and Blanchet and Schafer [Q. The second post- Newtonian correction, at 
order v 4 , was calculated in 1994 by Blanchet, Damour, Iyer, Will, and Wiseman |27||. Finally, the spin-orbit 
and spin-spin corrections were calculated in 1993 by Kidder, Will, and Wiseman ]2l|. 



We see from Eq. (31 ) that the post-Newtonian corrections bring a number of additional source parameters 
into the picture. The waveform no longer depends uniquely upon the chirp mass M,; it now depends upon the 
masses mi and 777.2 separately, and upon the spin-orbit and spin-spin parameters (which stay approximately 
constant as the system evolves toward coalescence). These new parameters must be included into /J when 
the signal is analyzed using matched filtering. 

What we have at this point is an expression for the waveform which is accurate to second post-Newtonian 
order. The question facing us is whether this waveform is sufficiently accurate to be used as measurement 
templates. In other words, are the systematic errors generated by these templates guaranteed to be smaller 
than the statistical errors? 

Evidently, to answer this question is difficult, since we do not have access to the exact waveform in order 
to make comparisons. In the next section we will consider a model problem for which the waveform can be 
calculated exactly, thereby enabling us to judge the accuracy of the post-Newtonian expansion. We will find, 
in Sec. 15, that the answer to this question is, most likely, no: the second post-Newtonian waveform is not 
sufficiently accurate. 



13 Waveform calculations: perturbation theory 

A different line of attack for solving Einstein's equations for a compact binary system is to assume that one 
of the bodies is very much less massive than the other p8| . We therefore demand 

fi/M < 1, (34) 

where \i is the reduced mass and M the total mass. In this limit \x is practically equal to the smaller mass 
777,1, and M is practically equal to the larger mass 777,2. In contrast with the post-Newtonian approach, we 
assume nothing about the size of the velocity v. 
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This approach takes advantage of the fact that when Eq. ( pi) ) is valid, the smaller mass creates only a 
small perturbation in the gravitational field of the larger mass, which can be taken to be the Schwarzschild 
field. The total gravitational field can therefore be written as 

9aP = 9a} + h ocf3, (35) 

where g^l represents the background Schwarzschild metric, and h a p the perturbation. 

Writing Einstein's equations for g a p and linearizing with respect to /i a /3, one finds that the perturbation 
must satisfy an inhomogeneous wave equation in the Schwarzschild spacetime. Schematically, 

W^h^ = SttT"' 3 , (36) 

where n a/3fJ - 1 ' is an appropriate curved-spacetime wave operator, and T a @ the stress-energy tensor associated 
with the orbiting mass. 

We will specifically assume that the central body is a Schwarzschild black hole of mass M . This as- 
sumption is made for simplicity, and removes the need to model the star's interior. As a matter of fact, 
the internal structure of the bodies is irrelevant, except during the last few orbital cycles before coalescence; 
this was discussed in Sec. 5. Taking advantage of this, we model the orbiting body as a point particle of 
mass /i, so that its stress-energy tensor is a Dirac distribution with support on the particle's world line. 
For simplicity, and also because it is physically well motivated (as explained in Sec. 5), we take the world 
line to be a circular geodesic of the Schwarzschild spacetime. With {t, r, 8, <fi} as the usual Schwarzschild 
coordinates, rg denotes the orbital radius, and £1 = dcj)/dt is the angular velocity. We have 

v = nr = (M/r ) 1/2 = (Mft) 1 / 3 , (37) 

and the gravitational-wave frequency / is given by 2irf = 20. We stress once more that in the perturbation 
approach, v is not required to be small. The only limitation on v comes from the fact that for < 6M, 
circular orbits are no longer stable; this implies v < 6 -1 / 2 ~ 0.4082. 

Black- hole perturbations are conveniently treated with the Teukolsky formalism |2j|, in which gravi- 
tational perturbations are represented by the complex- valued function ^4, a particular component of the 
perturbed Weyl tensor; the tensor h a p can be reconstructed from ^4. The equation satisfied by '5 4 admits 
a separation of the variables. When 'J 4 is expanded in spherical harmonics and decomposed into Fourier 
modes e~ luJt , one obtains an ordinary second-order differential equation — the Teukolsky equation — for 
the radial function Rwem{r). Here, £ and m are the usual spherical-harmonic indices. Schematically, and 
omitting the subscript ui£m, this equation takes the form 

VR{r) = T(r), (38) 

where T> is a second-order differential operator, and T(r) the source, constructed from the particle's stress- 
energy tensor p8| . 

Equation ( |38| ) can be integrated in the standard way by constructing a Green's function G(r,r') out 
of two linearly independent solutions to the homogeneous equation. These solutions, R<(r) and R > (r), 
respectively satisfy appropriate boundary conditions at the inner (r = 2M) and outer (r = 00) boundaries. 
Schematically, G(r,r') — R < (r < )R > (r > ), where r< (r>) is the lesser (greater) of r and r'. The solution to 
Eq. ( |38| ) can then be expressed as R(r) — J G(r, r')T(r') dr', and ^4 can be reconstructed by summing over 
all the modes. Finally, the gravitational waveform h and the luminosity dE/dt can be obtained from the 
asymptotic behaviour of ^4 when r — > 00. 

To integrate Eq. ( pq ) therefore reduces to solving the homogeneous Teukolsky equation, T>R(r) = 0, for 
the functions i?<(r) and R>(r). This, it turns out, is equivalent to integrating the Regge- Wheeler equation 
@ 

id 2 2 / 2M\ 



'£{£ + !) 6M~ 



>X uji (r) = 0, 



(39) 
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Figure 2: Various representations of (dE/dt)/(dE/dt)jy as a function of orbital velocity v. The solid curve 
represents the exact results, obtained numerically. The various broken curves represent the various post- 
Newtonian approximations, as explained in the text. 



for the functions -X"<(r) and X > (r); here, d/dr* = (1 — 2M/r)d/dr. This comes about because a solution 
to the homogeneous Teukolsky equation can easily be related to a solution to the Regge- Wheeler equation. 
The relation is known as the Chandrasekhar transformation fell. 



14 Luminosity from the perturbation approach 

The problem of calculating the gravitational waveform for a compact binary system with small mass ratio 
can therefore be reduced to the simple one of integrating Eq. (|39|). The relevant dimensionless parameter 
entering this equation is Mcu, where, for circular orbits, to = mfi [^8| . We therefore have, using Eq. j37|), 
Mto = mil 3 . 

For arbitrary values of Mui the Regge- Wheeler equation must be integrated numerically |32j . This must 
be done separately for each selected value of v in the interval 0.1 < v < 0.4 (approximately corresponding, 
for systems of a few solar masses, to the frequency interval 10 Hz < / < 1 000 Hz). H 

Once again we will focus on dE/dt, the gravitational- wave luminosity. In Fig. 2 we present a plot of 
(dE / dt) / {dE / dt) m as a function of v. This is the ratio of the luminosity as calculated exactly (numerically) 
using the perturbation approach, to the quadrupole-formula expression 




The numerical results are depicted as a solid curve. Apart from negligible numerical errors, these results are 
exact to all orders in v; the only assumption used in the calculation concerns the smallness of fi/M . The 
various broken curves will be described below. The figure shows that (dE / dt) / (dE / dt) m tends toward unity 
as v — ► 0, and stays within 15% of unity everywhere in the interval < v < 

It is now easy to judge the accuracy of the post-Newtonian expansion for dE/dt in the limit fi/M — > 0. 
One simply evaluates Eq. (J3l|) in this limit, with SO = SS = 0, and compare with the numerical results. This 
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post-Newtonian curve is labeled "4" in the figure, and we see that the numerical results are only imperfectly 
reproduced. The consequences of this will be discussed in the next section. 

One can in fact do better than this. The perturbation approach is not only suitable for exact, numerical 
computations. By combining it with a slow-motion approximation — putting w<lon top of fx/M <C 1 — it 
also becomes suitable for approximate, analytical computations. Methods for integrating the Regge- Wheeler 
equations analytically in the limit Mu> <C 1 were devised by various authors p4| p8| , |3^ |. The basic idea is 
to proceed by iterations: Eq. ( (39| ) with Mui — is solved in terms of spherical Bessel functions, and this 
zeroth order solution is used as input for the first iteration. 

Using these methods, Tagoshi and Sasaki |m| were able to calculate dE/dt analytically, accurately through 
fourth post-Newtonian order, two full orders beyond Eq. (|3l|). They derive the rather impressive expression 



dE 
~dt 
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(41) 



Notice the presence of lnw terms in this expansion, as well as that of the Euler number 7 ~ 0.5772. Notice 
also that the first four terms reproduce Eq. (|3l| ) in the limit fi/M ~ * 0, when SO = SS = 0. This, of course, 
is as it should be. 

The broken curves in Fig. 2 represent plots of Eq. (^) truncated to various orders in v. For example, 
the curve labeled "6" is a plot of Eq. (41) with all terms of order v 7 and v 8 discarded. We see from the 
figure that the post-Newtonian expansion converges poorly. (The suspicion, in fact, is that the series is only 
asymptotic.) Witness in particular the poor quality of the curve "5" with respect to "4", and compare also 
"8" to "7". 



15 Accuracy of the post-Newtonian expansion 

We have seen that in the /i/M -> limit, the post-Newtonian expansion for dE/dt converges poorly. A 
similar statement can be made about dE/df: the post-Newtonian expansion 

d_E(d_E\ f 3 81 675 6 _ 19845 \ 

df \df) N { * 8 16 128 + J' (42) 



where 



converges slowly to the exact result [p5| 



§) =~\,Mv-\ (43) 
df J N 3 



f-(f)(i-«^)(i-^. m 



' N 



Contrary to Eq. (fy]), the expansion (^2|) is actually known to converge (for all values of v in the interval 
< v< 6- 1 / 2 ). 

In this section we address the issue as to how much of an obstacle the poor convergence of the post- 
Newtonian expansion poses to the construction of accurate measurement templates J36|. We will answer 
this question for binary systems with small mass ratios, using the results described in the preceding section. 
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n 


PN 


EXACT 


4 


0.5796 


0.4958 


5 


0.4646 


0.5286 


G 


0.7553 


0.9454 


7 


0.7651 


0.9864 


8 


0.7568 


0.9695 



Table 1: Reduction in signal-to- noise ratio incurred when matched filtering with approximate, post- 
Newtonian templates. The first column lists the order n of the approximation, the second column lists 
lZ n as calculated using the post-Newtonian approximation for dE/df, and the third column lists TZ n as 
calculated using the exact expression for dE/df. 

Since there is no reason to believe that the convergence of the post-Newtonian expansion would be much 
improved for systems of comparable masses, our results should also apply, at least qualitatively, for such 
systems. 

One way to address this question is to ask, given a waveform constructed from the exact numerical 
results, how much signal-to-noise ratio is lost by matched filtering the exact signal with an approximate 
post-Newtonian template? 

From Sec. 8 we know that filtering with the exact signal h(t) would give SNR| max , the largest possible 
value of the signal-to-noise ratio. On the other hand, filtering with a post-Newtonian template h n , where n 
denotes the order in v to which the expansion is taken (for example, n = 4 represents a waveform accurate 
to second post-Newtonian order), gives the smaller value SNR| ac tuai- From the results of Sec. 7 and 8 we 
obtain 

n _ SNR| actual _ \(h\h n )\ 

" SNR| max y/(h\h)(h n \h n )' 

The detail of how to compute TZ n is presented in Ref . Q . 

The calculation described here can only be carried out for binary systems with small mass ratios, because 
only for these do we have access to the exact waveform h(t). Nevertheless, in the following we shall let 
n/M become large, without altering our expressions for the exact and post-Newtonian waveforms. This is 
done without justification, but reflects the adopted point of view that the quality of the post-Newtonian 
approximation should not be appreciably affected by the finite-mass corrections. 

In doing this extrapolation, some thought must be given as to the interpretation of the ratio ji/M. In 
the limit fi/M — * this can be taken to be both the ratio of the individual masses or the ratio of reduced 
mass to total mass. In the case of comparable masses, some choice must be made. We note that as fi/M 
is allowed to grow large, expressions ( ^0| ) and (|43|) for the Newtonian quantities (dE/dt)^ and (dE/df )n 
must be replaced by expressions (|2^) , which are formally identical. This shows that when extrapolating to 
the case of comparable masses, [i is to be interpreted as the reduced mass, and M as the total mass. 

We quote the results for lZ n corresponding to a system of two neutron stars, each of 1.4 M Q ; these are 
displayed in the second column of Table 1. We see that even at quite a high order in the post-Newtonian 
expansion, only three quarters of the signal-to-noise ratio is reproduced by the post-Newtonian template. 
This shows that the apparently small discrepancies between the exact and post-Newtonian results for dE/dt 
and dE/df provide a serious obstacle to the construction of accurate measurement templates. 

It is interesting to ask how much of the signal-to-noise ratio would be recovered if dE/df were kept exact 
instead of being expressed as a post-Newtonian expansion. The third column of the table displays the values 
of lZ n calculated in this way. We see that most of the signal-to-noise is recovered: TZ n can now be as large 
as 0.9864 instead of 0.7651. Why does the exact expression for dE/df give such better results? It can be 
established J|6| that this has to do with the following fact: While the exact expression for dE/df correctly 
goes to zero at v = 6 -1 / 2 (at the innermost circular orbit), its post-Newtonian analogue fails to do so. 
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[For example, the right-hand side of Eq. (f42|), truncated to order v 8 , goes to zero at v ~ 0.4236 > 6~ 1 / 2 , 
corresponding to a radius ro — 5.572M.] When corrected for this, the post-Newtonian template gives much 
better results. 

16 Conclusion 

We therefore see that the poor convergence of the post-Newtonian expansion is a serious obstacle to the 
construction of accurate measurement templates. Devising ways to extract the information contained in 
gravitational waves produced during the late inspiral of a compact binary system poses a great challenge to 
theoretical physicists. It is not clear that "simply" pushing to higher order in post-Newtonian theory will 
be enough to produce sufficiently accurate measurement templates (such that the systematic errors will be 
smaller than the statistical errors) . There may be a need for theorists to develop alternative ways of dealing 
with this problem. 

It should be stressed that the convergence problem does not arise when constructing detection templates. 
Indeed, as was discussed in Sec. 9, there is no particular need for these templates to be derived from the 
equations of general relativity. And since the detection templates need not involve any parameters of direct 
physical significance, the notion of systematic errors does not apply to them. The only requirement for 
constructing detection templates is that they should span the appropriate "signal space" , and that they 
should do so the most efficiently. 

The poor convergence of the post-Newtonian expansion is therefore not an obstacle for detecting gravita- 
tional waves from coalescing compact binaries. It is only an obstacle for extracting the information that the 
waves contain. To overcome this obstacle will undoubtedly be a theorist's challenge for years to come. This 
state of affairs is highly interesting from a historical point of view: Never before in the history of gravitational 
physics has experiment demanded such a high degree of sophistication on theoretical calculations. 
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